%function [P]=arma_spec(f,p,q,V);
%
%Computes the expected spectral distribution at frequencies, f, given
%autoregressive parametres, p; moving average parameters, q; and
%variance, v.  
%
%P. Huybers
%Harvard, 2021

function P=arma_spec(f,dt,p,q,v); 

if nargin<3,
    q=[];
end;
if nargin<4,
    v=1;
end;

pf=zeros(size(f));
for ct=1:length(p); 
    pf=pf+p(ct)*exp(-i*2*pi*f*ct*dt);
end;
qf=zeros(size(f));
for ct=1:length(q); 
    qf=qf+q(ct)*exp(-i*2*pi*f*ct*dt);
end;

%P=v*abs((1+qf).^2)./abs((1-pf).^2);
P=v*abs((1+qf)./(1-pf)).^2;
